function similarity_transforms()
% SIMILARITY_TRANSFORMS 相似变换演示
% 
% 展示相似变换在特征值计算中的作用和应用

fprintf('=== 相似变换演示 ===\n\n');

%% 1. 相似变换基础
fprintf('1. 相似变换基础\n');
fprintf('目标：理解相似变换的基本性质\n\n');

% 原始矩阵
A1 = [2, 1, 0; 1, 3, 1; 0, 1, 2];
fprintf('原始矩阵 A1 =\n');
disp(A1);

% 相似变换矩阵
P1 = [1, 2, 1; 0, 1, 1; 1, 0, 2];
fprintf('变换矩阵 P1 =\n');
disp(P1);

% 相似变换：B = P^(-1) * A * P
B1 = P1 \ A1 * P1;
fprintf('相似变换后的矩阵 B1 = P1^(-1) * A1 * P1 =\n');
disp(B1);

% 验证特征值不变性
eigenvals_A1 = sort(eig(A1));
eigenvals_B1 = sort(eig(B1));
fprintf('A1的特征值: [%.6f, %.6f, %.6f]\n', eigenvals_A1);
fprintf('B1的特征值: [%.6f, %.6f, %.6f]\n', eigenvals_B1);
fprintf('特征值差异: %.2e\n', norm(eigenvals_A1 - eigenvals_B1));

% 验证迹和行列式不变性
trace_A1 = trace(A1);
trace_B1 = trace(B1);
det_A1 = det(A1);
det_B1 = det(B1);
fprintf('迹: A1 = %.6f, B1 = %.6f, 差异 = %.2e\n', trace_A1, trace_B1, abs(trace_A1 - trace_B1));
fprintf('行列式: A1 = %.6f, B1 = %.6f, 差异 = %.2e\n', det_A1, det_B1, abs(det_A1 - det_B1));

%% 2. 对角化变换
fprintf('\n2. 对角化变换\n');
fprintf('目标：通过相似变换将矩阵对角化\n\n');

A2 = [4, 2; 1, 3];
fprintf('矩阵 A2 =\n');
disp(A2);

% 计算特征向量矩阵
[V2, D2] = eig(A2);
fprintf('特征向量矩阵 V2 =\n');
disp(V2);
fprintf('对角矩阵 D2 =\n');
disp(D2);

% 验证对角化：A = V * D * V^(-1)
A2_reconstructed = V2 * D2 * V2^(-1);
fprintf('重构误差: %.2e\n', norm(A2 - A2_reconstructed, 'fro'));

% 相似变换到对角形式
D2_check = V2^(-1) * A2 * V2;
fprintf('通过相似变换得到的对角矩阵:\n');
disp(D2_check);
fprintf('与eig结果的差异: %.2e\n', norm(D2 - D2_check, 'fro'));

%% 3. Householder变换
fprintf('\n3. Householder变换\n');
fprintf('目标：使用Householder变换进行相似变换\n\n');

A3 = [4, 3, 2, 1; 3, 3, 2, 1; 2, 2, 2, 1; 1, 1, 1, 1];
fprintf('对称矩阵 A3 =\n');
disp(A3);

% 第一步Householder变换：消除第一列的下三角元素
x = A3(2:end, 1);  % 需要消除的向量
fprintf('需要消除的向量 x = [%.2f, %.2f, %.2f]''\n', x);

% 计算Householder向量
alpha = -sign(x(1)) * norm(x);
u = x;
u(1) = u(1) - alpha;
u = u / norm(u);

fprintf('Householder向量 u = [%.4f, %.4f, %.4f]''\n', u);

% 构造Householder矩阵 H = I - 2*u*u'
n = length(x);
H_small = eye(n) - 2 * (u * u');
H3 = eye(size(A3));
H3(2:end, 2:end) = H_small;

fprintf('Householder矩阵 H3 =\n');
disp(H3);

% 应用相似变换
A3_transformed = H3' * A3 * H3;
fprintf('变换后的矩阵 A3_transformed =\n');
disp(A3_transformed);

% 检查是否成功消除了目标元素
fprintf('第一列下三角元素: [%.2e, %.2e]\n', A3_transformed(3,1), A3_transformed(4,1));

% 验证对称性保持
symmetry_error = norm(A3_transformed - A3_transformed', 'fro');
fprintf('对称性误差: %.2e\n', symmetry_error);

%% 4. Givens旋转变换
fprintf('\n4. Givens旋转变换\n');
fprintf('目标：使用Givens旋转进行相似变换\n\n');

A4 = [3, 2, 1; 2, 4, 2; 1, 2, 3];
fprintf('对称矩阵 A4 =\n');
disp(A4);

% 消除A4(3,1)元素
a = A4(1, 1);
b = A4(3, 1);
fprintf('目标：消除元素 A4(3,1) = %.4f\n', b);

% 计算Givens旋转参数
if abs(b) < 1e-12
    c = 1; s = 0;
else
    if abs(b) > abs(a)
        t = -a/b;
        s = 1/sqrt(1 + t^2);
        c = s * t;
    else
        t = -b/a;
        c = 1/sqrt(1 + t^2);
        s = c * t;
    end
end

fprintf('Givens参数: c = %.6f, s = %.6f\n', c, s);

% 构造Givens矩阵
G4 = eye(3);
G4(1, 1) = c; G4(1, 3) = -s;
G4(3, 1) = s; G4(3, 3) = c;

fprintf('Givens矩阵 G4 =\n');
disp(G4);

% 应用相似变换
A4_rotated = G4' * A4 * G4;
fprintf('旋转后的矩阵 A4_rotated =\n');
disp(A4_rotated);

% 检查目标元素是否被消除
fprintf('目标元素 A4_rotated(3,1) = %.2e\n', A4_rotated(3,1));

% 验证正交性
orthogonality_G4 = norm(G4' * G4 - eye(3), 'fro');
fprintf('G4的正交性误差: %.2e\n', orthogonality_G4);

%% 5. 三对角化变换
fprintf('\n5. 三对角化变换\n');
fprintf('目标：将对称矩阵变换为三对角形式\n\n');

A5 = [4, 1, 2, 1; 1, 3, 1, 2; 2, 1, 5, 1; 1, 2, 1, 4];
fprintf('对称矩阵 A5 =\n');
disp(A5);

% 使用MATLAB内置函数进行三对角化
[Q5, T5] = hess(A5);
fprintf('正交变换矩阵 Q5 =\n');
disp(Q5);
fprintf('三对角矩阵 T5 =\n');
disp(T5);

% 验证变换
A5_reconstructed = Q5 * T5 * Q5';
fprintf('重构误差: %.2e\n', norm(A5 - A5_reconstructed, 'fro'));

% 检查三对角结构
upper_diag = diag(T5, 1);
main_diag = diag(T5);
lower_diag = diag(T5, -1);
fprintf('主对角线: [%.4f, %.4f, %.4f, %.4f]\n', main_diag);
fprintf('上对角线: [%.4f, %.4f, %.4f]\n', upper_diag);
fprintf('下对角线: [%.4f, %.4f, %.4f]\n', lower_diag);

% 检查其他元素是否为零
other_elements = T5;
other_elements = other_elements - diag(main_diag) - diag(upper_diag, 1) - diag(lower_diag, -1);
fprintf('非三对角元素的范数: %.2e\n', norm(other_elements, 'fro'));

%% 6. Hessenberg变换
fprintf('\n6. Hessenberg变换\n');
fprintf('目标：将一般矩阵变换为Hessenberg形式\n\n');

A6 = [1, 2, 3, 4; 5, 6, 7, 8; 9, 10, 11, 12; 13, 14, 15, 16];
fprintf('一般矩阵 A6 =\n');
disp(A6);

% Hessenberg变换
[Q6, H6] = hess(A6);
fprintf('正交变换矩阵 Q6 =\n');
disp(Q6);
fprintf('Hessenberg矩阵 H6 =\n');
disp(H6);

% 验证变换
A6_reconstructed = Q6 * H6 * Q6';
fprintf('重构误差: %.2e\n', norm(A6 - A6_reconstructed, 'fro'));

% 检查Hessenberg结构（下三角部分除了次对角线应为零）
lower_part = tril(H6, -2);
fprintf('次对角线以下元素的范数: %.2e\n', norm(lower_part, 'fro'));

%% 7. 平衡变换
fprintf('\n7. 平衡变换\n');
fprintf('目标：改善矩阵的数值性质\n\n');

% 构造不平衡矩阵
A7 = [1e-3, 1; 1e3, 1];
fprintf('不平衡矩阵 A7 =\n');
disp(A7);

% 计算条件数
cond_A7 = cond(A7);
fprintf('原始条件数: %.2e\n', cond_A7);

% 平衡变换
[T7, B7] = balance(A7);
fprintf('平衡变换矩阵 T7 =\n');
disp(T7);
fprintf('平衡后的矩阵 B7 =\n');
disp(B7);

% 计算平衡后的条件数
cond_B7 = cond(B7);
fprintf('平衡后条件数: %.2e\n', cond_B7);
fprintf('条件数改善比: %.2f\n', cond_A7 / cond_B7);

% 验证相似性
eigenvals_A7 = sort(eig(A7));
eigenvals_B7 = sort(eig(B7));
fprintf('特征值保持不变，差异: %.2e\n', norm(eigenvals_A7 - eigenvals_B7));

%% 8. 相似变换的应用：矩阵幂的计算
fprintf('\n8. 相似变换的应用：矩阵幂的计算\n');
fprintf('目标：通过对角化高效计算矩阵幂\n\n');

A8 = [2, 1; 1, 2];
k = 10;  % 计算A^10
fprintf('矩阵 A8 =\n');
disp(A8);
fprintf('计算 A8^%d\n', k);

% 方法1：直接计算（效率低）
tic;
A8_power_direct = A8^k;
time_direct = toc;

% 方法2：通过对角化计算
tic;
[V8, D8] = eig(A8);
D8_power = D8^k;  % 对角矩阵的幂很容易计算
A8_power_eigen = V8 * D8_power * V8^(-1);
time_eigen = toc;

fprintf('直接计算结果:\n');
disp(A8_power_direct);
fprintf('对角化计算结果:\n');
disp(A8_power_eigen);
fprintf('结果差异: %.2e\n', norm(A8_power_direct - A8_power_eigen, 'fro'));

fprintf('计算时间比较:\n');
fprintf('直接方法: %.6f 秒\n', time_direct);
fprintf('对角化方法: %.6f 秒\n', time_eigen);

% 对于大幂次，对角化方法的优势更明显
k_large = 100;
fprintf('\n计算 A8^%d:\n', k_large);

tic;
D8_power_large = diag(diag(D8).^k_large);
A8_power_large = V8 * D8_power_large * V8^(-1);
time_large = toc;

fprintf('对角化方法计算 A^%d 用时: %.6f 秒\n', k_large, time_large);

%% 总结
fprintf('\n=== 相似变换总结 ===\n');
fprintf('1. 相似变换保持特征值、迹、行列式不变\n');
fprintf('2. 对角化是最重要的相似变换\n');
fprintf('3. Householder变换适用于对称矩阵的三对角化\n');
fprintf('4. Givens旋转可以精确控制零化元素\n');
fprintf('5. 三对角化简化特征值计算\n');
fprintf('6. Hessenberg形式适用于一般矩阵\n');
fprintf('7. 平衡变换改善数值稳定性\n');
fprintf('8. 相似变换在矩阵函数计算中很有用\n\n');

end